Scaling Invariance and the Iterative 
Transformation Method for a Class of 
Parabolic Moving Boundary Problems 

Riccardo Fazio 
Department of Mathematics and Computer Science, 
University of Messina, Viale F. Stagno D'Alcontres, 31 
98166 Messina, Italy 
E-mail: rfazio@unime.it 



Home-page: http://mat52 1 .unime.it/fazio 



November 29, 2012 

Abstract 

In this paper we apply a scaling invariance analysis to reduce a class of 
parabolic moving boundary problems to free boundary problems governed 
by ordinary differential equations. As well known free boundary problems 
are always non-linear and, consequently, their numerical solution is often 
obtained iteratively. Among the numerical methods, developed for the nu- 
merical solution of this kind of problems, we focus on the iterative trans- 
formation method that has been defined within scaling invariance theory. 
Then, as illustrative examples, we solve two problems of interest in the ap- 
plications. The obtained numerical results are found in good agreement with 
exact or approximate ones. 

Key Words. Scaling invariance, numerical transformation method, Stefan's prob- 
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1 Introduction. 



The main contribution of this paper is the development of a complete scaling in- 
variance analysis for a class of parabolic moving boundary problems and apply an 
extended scaling transformation to solve the reduced free boundary problems by 
an initial value solver. In fact, the scaling invariance analysis allows us to reduce 
the original class of problems to a class of free boundary problems governed by 
ordinary differential equations. As pointed out first by Landau G3ll . free or mov- 
ing boundary problems are always non-linear, and therefore they are often solved 
numerically. In fact, a superposition principle for the solution of these problems 
cannot be valid because a variation of the auxiliary data, initial or free boundary 
conditions, produces a change in the free or moving boundary that in turn changes 
the domain of existence of the solution. 

Most of the existing references on this subject, such as the one in the book by 
Crank Q, or the one by Tarzia GTTl . are essentially devoted to the solution of the 
well known Stefan problem which is governed by the linear heat equation. For the 
numerical solution of such a problem several different approaches have been con- 
ceived over the years. Among the most famous ones there are the front-tracking, 
the front-fixing, and the domain-fixing methods (see [7, pp. 217-281]), as well as 
other finite-difference or finite element approaches (see, for instance, Meek and 
Norbury [|24ll . Bonnerot and Jamet [61 or Asaithambi [0), or moving grid, level 
set, or phase field methods (see the review by Javierre [22J). Unfortunately, the 
proposed methods are introduced in the case of linear parabolic partial differential 
equations and they are not easily extended to non-linear parabolic cases belonging 
also to the class of problem (Q~|). In this paper, we propose a method to overcome 
these difficulties, provided that the problem is invariant with respect to a scal- 
ing group. To this end, we use the similarity approach, described in full details 
by Dresner [8], within Lie's group invariance theory (see Bluman and Cole flU, 
Dresner [9], Barenblatt Q, or Bluman and Kumei [0). As far as the performance 
of different methods is concerned, the introductory remark in a survey paper by 
Fox 11201 is pertinent: "Problems of the same general nature can differ enough in 
detail to make a good method for one problem less satisfactory and even mediocre 
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for another almost similar problem". This point of view justifies the development 
of so many different numerical methods. 

The application of scaling invariance to applied mathematics and numerical 
analysis has been a fruitful research field for more than a century, see Fazio IfToH . 
As far as numerical applications are concerned, the first numerical transformation 
method (TM) is due to Topfer Il28ll . He solved non-iteratively, using a transforma- 
tion of variables, the Blasius problem of boundary layer theory. Several problems 
in boundary-layer theory lack this kind of invariance and cannot be solved by non- 
iterative (I)TMs, see [25, Chapters 7-9]. To overcome this drawback it is possible 
to define an iterative extension of the Topfer' s algorithm lfT8l [T2l [TTl [TBI . Nu- 
merical solution of free boundary problems by means of TMs, developed within 
scaling invariance theory, are considered in |[T8l [TOl [T5ll . Here, we explain in 
full details the definition of the ITM. Moreover, in order to show the validity of 
the proposed approach we solve two relevant problems: the single phase Stefan's 
problem ll26l [71 EJJ, and a problem describing the spreading of a viscous fluid 
above a smooth horizontal surface |[24l . For the former we compare the results 
obtained by the proposed ITM with those given by an asymptotic analysis. As far 
as the second problem is concerned, the exact similarity solution is available, and 
therefore we are able to present a direct test for the obtained numerical results. 



2 Scaling invariance 

We consider the following class of moving boundary problems of the parabolic 
type 



du 
~dt 



— \ 

dx 




on t > 0, < x < x w (t) , 



u(x,0) = , 




u(0,t) =At 



a 



(1) 
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dx 



(Xw(t),t) 



t,x w (t),^-(t),u(x w (t),t) 



where n, A, a, B, and j8 are constants, x and t represent time and space respec- 
tively, u(x,t) is the field variable, /?(-,•?■)> an d <?(•,-,■,•) are given functions of 
their arguments, and x w (t) is the unknown moving boundary. 

As mentioned before, the problem (GQ) is non-linear because x w (t) depends on 
the initial and boundary data so that a superposition principle cannot be valid (that 
was pointed out by Landau 11231 ). As a consequence, obtaining analytical solutions 
for problems belonging to the class © is a difficult task (see [7, pp. 101-139]). 

The following scaling group 



(2) 



where A is the (positive) group parameter, leaves the problem (QQ) invariant pro- 
vided that 
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7- 



nOC + 1 

p(-,v 



aut a = - ^ + 1 . with 6 / - - 1 
2 — n — np n 



t a p(x w {t)t- l lv, 



dt 



!H( f ) f (r-i)/y 



(3) 



q(-,., ■, ■) = t^-^Q (x^t-Vr, ^(t)t^l^u{x w {t),t)r a 

As a consequence, we can introduce the similarity variables as follows: 



7] =Xt 



-1/y 



Tj w — X w (t )? 



-1/y 



U(rj)=t a u(x,t) . 



(4) 



By using ©, we see that the model problem (OQ) reduces to 



d 2 U 
drj : 



dU 



dr\ 

U(0) =A 



dU 



+ nU~ I — — +-t]U~ n - aU 



rl-n 
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dr\ 



aut 4^(0) =5 



(5) 



£/(t? w ) =P(t?w,t?w) , j^iVw) =Q(t? w ,tj w ,£/(t? w )) , 

where 7] w is the unknown free boundary for the ordinary differential problem ©. 
We notice that, for any functional form of P(-,-) and Q(-,-,-), the free bound- 
ary conditions in © depend only on r\ w . Figure Q] shows the map given by the 
similarity variables ©. 
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3 The numerical method 



Here we consider the class of free boundary problems 

d 2 w „ / dw\ 



dz 2 f \ z ' w 'dz)' (6) 
*(w(0),^(0))=C, w(s)=j(s) 1 j z {s)=t{s), 

where w and z are the field and independent variable, respectively, /(•, •, •), g(-, •), 
j(-), and £(■) are given functions of their variables, C is a given constant and s 
represents the unknown free boundary. Of course, it is a simple matter to verify 
that © belongs to © for appropriate choices of /, g, j, and t. A free boundary 
problem belonging to © can be solved by the ITM defined by the following steps. 

-) First we introduce an extended problem, namely: 

= / ( h- S /°Z,h-V°wM d - 1)/a ^-) , 

dz 1 V dz J 

h l l°g (h-V°w(0)M S - 1)/a ^(0)) = C, (7) 

w(s) = h 1 / a j(h- s / (7 s) , — (h- 8 l a s) = h^- s y (7 £(h- s / (y s) , 

dz 

where C 7^ 0, and h is a parameter. A constructive characterization of CO), 
within similarity analysis, is given in |[T4ll . Let us remark that the free 
boundary problem © is recovered from the extended problem © by set- 
ting h = 1 . Moreover, the extended problem © is partially invariant with 
respect to the extended scaling group 

z* = (0 8 z , s* = co s s , w* = cow , h* = co a h , (8) 

where 00 is the (positive) group parameter, while 8 and o are constants re- 
lated to the particular problem under study. We notice that the governing 
equations and the two free boundary conditions are invariant, but the condi- 
tion at z = is not invariant. 

-) Given the values of 8, and h*, we fix a value of s* greater than zero, and 
integrate ©, written in the starred variables, inwards in [0,s*] to compute 
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approximate values of w*(0) and m order to get 

h *\/a g ( h *-l/a w *( ^ h *(8-l)/o^(fy\ 
(0 = ^ z (9) 

and, by using the scaling invariance, the corresponding value ofh = (0~ a h*. 

-) We get a solution of the original free boundary problem © when we find a 
value of h* that transforms to h = 1. This is equivalent to find a zero of the 
so called transformation function 

T(h* ) = co^h* - 1 . (10) 

To this end we can apply a root finder or a bracketing method. The values 
of interest are defined by the scaling relations: 

*=arV, w(0) = io-V(0) > ^(0) = CO 5 ' 1 ^-(0). (11) 

Within the iteration we define the sequences h* and for j = 0,1,2, 

If r(/ip tends to zero as y goes to infinity, then Sj goes to the correct free 
boundary value s in the same limit. 

For both examples reported in the next section we used, as initial value solver, 
the classical fourth order Runge-Kutta scheme and the secant method as root 
finder. The convergence criterion for the secant method, was given by 

\r(h*j) | < Tol , and \sj - | < Tol , (12) 

where Tol = ID — 06. Here and in the following, the D notation indicate a double 
precision arithmetic. 



4 Two applications 

As a first example, we consider the classical parabolic moving boundary prob- 
lem: the celebrated Stefan's problem Q pp. 2-4, p. 9]. This problem, in non- 
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dimensional variables, can be written as follows 

du d 2 u ^ 

-=— = on ?>0, 0<x<x w (t), 

ot ax L 

h(x,0) = 0, Xw (0) = 0, 

m(0,0 = 1, (13) 

u(x w (t),t) = , 
du . . . . „dx w . . 

_(,„(,),<) = -^w- 

The last boundary condition, involving the free boundary velocity, is known as the 

Stefan's equation, and 1/5 is the Stefan's number. The scaling invariance analysis 

of section 2 is valid for £□]) by setting n = 0, A = 1, a = 0, y = 2, P(-, •) = 0, 
dx 

and Q(',-,-) = ~^^ 2 ~^' Therefore, by introducing the similarity variables ©, 
with a = and y = 2, the problem (TT3T) reduces to 



J 2 [/ 1 rftf _ 

f/(0) = 1 , (14) 

Newmann's solution can be easily expressed by using the special error function 
erf(-), see Q, 

U(l) = 1 ~ erf(T?/2) (15) 
The free boundary, T] w , is the positive real root of the equation 

^ 1 / 2 577 w exp(7] vl 2 /4)erf(7] w /2) -2 = 0. (16) 

This equation is transcendental and its solution allows us to obtain the exact mov- 
ing boundary solution. The solution of (fT6l) is unique ,and this implies the unique- 
ness of the similarity solution. As a consequence, the Stefan's problem admits 
only one solution. Approximate values of r} w , for different values of S, can be 
found by a standard error function asymptotic expansion, see II2TTI . 

In the following we want to find approximate values of r\ w by using the ITM. 
First we have to introduce an extended problem 

d 2 U h 1 / 2 dU__ 
^7 2 " + ^~ r? ^7"" ' 



[7(0) = 1 , (17) 

dU /i 3/4 
U(7] w )=0, j^(Vw) = — 2~Sti w , 

and the extended scaling group 

rj* = co 1 ri, n w * = (o~ l ri w , U* = coU, h* = a 4 h . (18) 
So that, by using the scaling invariance, it follows that 

ITT IT T5f: 

G) = U*(0), h = (0- 4 h*, — (0) = AT 2 -— (0) , 7] w = (07] w *. (19) 

Let us remark that, due to the h 1 / 2 term in (fTTT) . we are allowed to consider positive 
values of h* only. In table Q] we list numerical results obtained for several values 
of S. Comparing the data reported in the last two columns of table \T\ we see that 
there is a good agreement between the results obtained by the present approach 
and those obtained by the asymptotic one. For the sake of brevity we omitted to 
list the intermediate data for the reported iterations. The secant method always 
verified the convergence criterion (fT2l) in few iterations. Figure [2] shows a sample 
numerical solution. 

As a second example, we consider a problem describing the spreading of a vis- 
cous fluid, such as a teacle, under the action of gravity above a smooth horizontal 
surface, see [[24]. Let us introduce the mathematical model 



du d 
dt dx 



3 du 
dx 



on t > , < x < x w (t) 



m(jc,0)=0, x w (0)=0, 

|^(0,0 = 0, (20) 

u(x w (t),t) =Ht~ l l 5 , 

^-{x w {t),t) = Lx w ~ x ^^{t)u{x w {t),t)~ 3 ' , 
ox dt 

where H and L are given constants. The scaling invariance analysis of section 2 is 

valid for ([20]) by setting n = 3, B = 0, j8 = -1/5, y= 5, P(-, ■) = H, and £>(•, •, •) = 
dx 

Lt 2 / 5 x w ~ l — -(t)u(x w (t), t)~ 3 . Therefore, by using the similarity variables ©, the 
dt 
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free boundary problem for (1201) is given by 

d 2 U „ , fdU\ 2 1 ,d£/ 1 2 „ 

^(0) = 0, (21) 

f(l*)=», ^(U = ^J- 
We notice that the boundary condition at r\ = is homogeneous, that is C = in 
©, and, as a consequence, in order to apply the ITM, we have to introduce a new 
dependent variable, namely 

V(r 1 ) = U(r 1 ) + r 1 . (22) 
Hence, the free boundary (|2TT) becomes 
d 2 V -ifdV \ 2 1 ^sfdV \ 1 2 



£(0)-l, (23) 



In order to apply the ITM, as a first step, we have to introduce an extended problem 

1/2 N 



dr\- 



. + 3C v-^)-g-*'/ 



+y H(V - * 1/2 »1)- 3 -* 1/2 ) + f (V - A'' 2 ir 2 = , 

£(0)=1, (24) 

V(l)„)=M?+A 1 <V, Sy(''»' = ' , ' /2 (^|3+ 1 ) ' 
and the extended scaling group 

7? * = £0 1 /2 77) Th,* = fl) 1 / 2 ^ , V* = COV, h* = ah. (25) 

Applying the scaling invariance properties we get 

dU 
dri 



£0=1-— (0)) , h = a)- l h\ £7(0) = co"V*(0) , 7] w = co 1/2 ?7 W * . 



(26) 
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Now, let us compare our numerical results with an exact similarity solution avail- 
able in literature for specific values of the parameters. In fact, in the particular 
case where H = 1/2 and L = — 1/2 the exact solution 



U{r\) 



3 ' 5 2 .2 



10 Vl2 +7? - 



1 1/3 

T? w = 1 , (27) 



was quoted in [|24l . In this case, from (1271) we get U (0) 0.751847, that is correct 
to six decimal places. Table[2]lists some numerical results obtained by the present 
approach for this particular choice of the parameter values. At it is easily seen, 
the reported values for U (0) and r\ w are in a very good agreement with the exact 
solution (1271) . Moreover, from the same table, we realize that comparable results 
must be found with different choices of T] w * and that, as far as the use of the secant 
method is concerned, we do not need to bracket the root of the transformation 
function. Figure [3] shows the obtained numerical solution. 

As far as the model ([201) is concerned, the fluid flux at x = is given by 

u(0,t) n ^(0,t)=Bt^ n+l ^[U(0)} n ^-(0) . (28) 

However, the problem considered in 11241 and studied herein for a viscous fluid 
has zero flux at x = because 5 = 0. On the other hand, the value of U (0) is of 
interest because it defines the fluid height at x = according to 

u(0,t) = t P U(0) , (29) 

where j8 = — 1/5 for the considered fluid. 



5 Concluding remarks 

This paper makes two main contributions: the scaling invariance analysis for the 
class of parabolic moving boundary problems (OQ) allows us to characterize those 
problems that can be reduced by similarity variables to free boundary problems 
governed by ordinary differential equations, and the definition of the ITM for the 
numerical solution of this second type of problems. From a numerical viewpoint, 
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it is simpler to solve a free boundary problem governed by an ordinary differen- 
tial equation than to face a moving boundary problem governed by a parabolic 
partial differential equation. The development of the ITM was motivated by the 
limitations of non-ITMs that sometimes can be applied due to the scaling invari- 
ance properties of the considered model, see Il25l Chapters 7-9]. The ITM has 
found application to the numerical solution of parabolic problems, either moving 
boundary problems or problems defined on infinite domains. In particular, in [fTTl 
the ITM is used to solve the sequence of free boundary problems obtained by a 
semi-discretization of ID parabolic moving boundary problems, and in |fT9l a free 
boundary formulation for the reduced similarity models is used in order to propose 
a moving boundary formulation for parabolic problems on unbounded domains. 

The present approach can be used for solving problems with different values 
of the parameters involving in (OQ) as well as different functional forms of 

P ^ w (0,^(o)=^(^(^- 1/r ,^(^ (r - 1)/r ) 

q^t,x w (t),^(t),u(x w (t),t)j = 

= ? («r-i)/r e ^ w(;)r i/r^ w? (r-i)/y M(Xw(f)jf)r 

Two problems, that have been defined in the applied sciences, were studied within 
the proposed framework. The numerical results reported in the previous section 
clearly show the correctness and reliability of our approach. 
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Figure 1: The similarity map. Note that the similarity lines x = and x w (t) are 
mapped to r\ = and rj w , respectively. 




Figure 2: Numerical solution for Stefan's problem with 5=1. We used the results 
reported in table [Hand the Runge-Kutta method with one hundred of steps. 
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Figure 3: Numerical solution for the problem (12TT) with L = 1/2 and H = —1/2. 
We notice that the exact similarity solution verifies the free boundary conditions: 
U (r/ w ) = 1/2 and ~j^( r \w) = —4/5. We used the results reported in the last line 
of tableland the Runge-Kutta method with one hundred of steps. 
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s j 



ITM 

*5 > 



Asymptotic 11211 

Tj w 7] w 



0.1 





600. 










1 


700. 










10 


639.263216 


-0.610425 


2.514145 


2.513961 


0.5 





100. 










1 


150. 










9 


105.180667 


-0.760017 


1.601231 


1.60H87 


1.0 





30. 










1 


40. 










8 


37.843777 


-0.910875 


1.240134 


1.240161 


5.0 





3. 










1 


2. 










7 


2.256999 


-l. 683000 


0.612848 


0.612864 


10.0 





1. 










1 


0.5 










8 


0.599873 


-2.309323 


0.440033 


0.440000 


50.0 





l.D-03 










1 


l.D-02 










11 


2.53D-02 


-5.03323 


0.199338 


0. 199499 



Table 1: Iterations obtained with rj w * = 0.5 and a step size Arj* = - ID -03. 
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j 




r(*5) 


1/(0) 




0.5 





0.5 


0.177999 








1 


0.1 


-0.198207 








7 


0.250158 


-3.14D-08 


0.751803 


0.996840 


1.0 





0.5 


-0.152895 








1 


0.1 


-0.349655 








7 


1.00032 


-1.78D-08 


0.751825 


0.999842 



Table 2: In the iterations we used a step size Ati* = — 5D — 04. 
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